Kaplan-Meier estimation.
Analyser og visualisering i R.
Introduktion til registerdata-sættet og hvordan i går til det resten af kurset.
Øvelser.
R-pakker relevant for forløbsanalyseanalyse.R.R.Kaplan-Meier overlevelsesanalyse i R (univariat og bivariat, inkl. test af forskelle).Overlevelse som funktion af tiden udtrykkes som\[ S(t) = Pr(T > t)\]
Hazard-raten udtrykkes som\[ h(t) \]
og
\[ h(t) \Rightarrow 0 \: \text{to} \: \infty \] da det er en rate og ikke en sandsynlighed.
Median leveltid udtrykkes som\[ S(t) = 0,5 \]
sandsynligheden for at tiden fra starten af forløbet indtil hændelsen er større end \(t\).
\[ P \left ( \left \{ T > t \right \} \right ) \]
Denne sandsynlighed kommer til udtryk som (integralet af) en sandsynlighedsfunktion.
\[ S(t) = P(\{T>t\}) = \int_{t}^{\infty} f(t)dt \]
sandsynligheden for at en begivenhed finder sted i et defineret tidsinterval, betinget af overlevelse op til det tidspunkt. Hazard-raten udtrykker dermed “styrken” af en begivenheds sandsynlighed, udtrykt som:
\[ h(t) = \frac{f(t)}{\hat{S}(t)} \]
Hvor \(f(t)\) er sandsynlighedsfunktion.
den tidligste tid hvor halvdelen af populationen/sample oplever begivenheden. Vi er altså interesseret i tiden \(t\) (den afhængige variabel i overlevelsesanalysen) når:
\[ \hat{S}(t) = 0,5\]
\(S(t)\) er et mål for overlevelse som en funktion af tiden. Med andre ord, hvor længe man bliver i en pågældende position, bestemt som overlevelsessandsynligheden i alle de foregående tider ganget sammen. Fortæller noget om den samlede udvikling.
\(h(t)\) fortæller noget om forekomsten af hændelser — forudsat personen har “overlevet” frem til tidsintervallet — samt hvordan det udvikler sig over tid. Fortæller noget om udviklingen i de enkelte tidsintervaller.
\(S(t) = 0,5\) fortæller noget om, hvornår halvdelen af sample/population ikke længere er under risiko for at opleve hændelsen. Fortæller noget om et enkelt punkt på kurven.
Med eksemplet fra øvelsen:
KM er en non-parametric metode til at estimere risikoen (sandsynligheden) for et udfald efter et givent tidspunkt baseret på de observerede “overlevelses”-tider (dur) i population (df) — i.e., der er ingen a priori antagelser om den underliggende fordeling af overlevelsestider):
\[ \hat{S}(t) = \prod_{i: t_i < t} \left(1 - \frac{d_j}{n_j}\right) \]
hvor, \(\hat{S}(t)\), overlevelsesfuntionen, er den estimerede sandsynlighed for overlevelse ved tid \(t\). \(t_{j}\) er observede tidspunkter for begivenheden (eng.: failure times). \(d_{j}\) er antallet af begivenheder ved tid \(t_{j}\). \(n_{j}\) er antallet af individer der ikke har oplevet begivenhed eller er blevet censureret før tid \(t_{j}\).
Med andre ord, vi multiplicere overlevelsessandsynlighederne over alle diskrete tidsintervaller (her, uger) hvor flyt finder sted. For hvert diskrete tidspunkt, \(t_{j}\), estimeres overlevelsessandsynligheden som produktet af overlevelsessandsynlighederne op til det tidspunkt, justeret for antallet af observerede events, \(d_{j}\), og antallet af personer i risiko for flyt, \(n_{j}\), op til det tidspunkt.
Det primære “bidrag” med denne metode er at al information frem til censurering anvendes, fremfor at droppe de individer/observationer, der er censureret på et senere tidspunkt; det står delvist i kontrast til typiske tilgange i klassisk regressionsanalyse.
Altså,
\(d_{j}\) er antallet af hændelser i det diskrete tidsinterval.
\(n_{j}\) er antallet af person i risikosættet. Dvs. antallet af personer, der ikke har oplevet hændelsen eller blevet censureret op til tid \(j\). Det vil sige at \(S_{t}\) er udledt af risikosættet: ofte refereret til som \(R(t)\).
\(1 - \frac{d_j}{n_j}\) er sandsynlighed for overlevelse, betinget af overlevelse frem til \(j\).
Med \(\prod_{i: t_i < t}\) ganger vi overlevelsessandsynligheder frem til det pågældende tidspunkt, \(j\).
For eksempel:
| Tid | Risikosæt \(n_{i}\) |
Hændelser \(d_{j}\) |
\(n_{j}-d_{j}\) | \(\frac{n_{j}-d_{j}}{n_{j}}\) | \(\hat{S}(t_{j})\) |
|---|---|---|---|---|---|
| 0 | 77 | 0 | 0 | 1 | \(\hat{S}(t_{j})=1\) |
| 1 | 77 | 4 | 73 | \(1 \times \frac{73}{77}\) | \(\hat{S}(t_{j})=0,95\) |
| 2 | 73 | 10 | 63 | \(0,95 \times \frac{63}{73}\) | \(\hat{S}(t_{j})=0,82\) |
| 3 | 73 | 63 | 14 | \(0,82 \times \frac{49}{63}\) | \(\hat{S}(t_{j})=0,64\) |
| 4 | 73 | 49 | 8 | \(0,64 \times \frac{41}{49}\) | \(\hat{S}(t_{j})=0,53\) |
| 5 | 73 | 41 | 6 | \(0,53 \times \frac{35}{41}\) | \(\hat{S}(t_{j})=0,45\) |
| \(t_{j}\) | … | … | … | … | \(\hat{S}(t_{j})=0\) |
\[ \hat{S}(4) = \frac{73}{77} \times \frac{63}{73} \times \frac{49}{63} \times \frac{41}{49} = 0,53 \]
\(\hat{S}(t)\) vil altid være stigende eller konstant. Aldrig faldende.
Alternativt kan KM estimeres som:
\[ \hat{S} (t_\left({i}\right)) = \hat{S} (t_\left({i-1}\right)) \times P \left ( T > t | T \geq t \right) \]
Hvilket udtrykker sandsynligheden for overlevelse indtil tidspunkt \(t_\left({j-1}\right)\) gange den betingede sandsynlighed for overlevelse for til \(t_\left({j}\right)\). For eksempel, ved \(t=4\):
\(\hat{S} (t_\left({4-1}\right)) = \frac{73}{77} \times \frac{63}{73} \times \frac{49}{63}\)
og
\(P \left ( T > 4 | T \geq 4 \right) = \frac{35}{41}\)
Et af de typiske formål med en KM estimation, er at sammeligne to overlevelseskurver. Vi gør dette ved at tilføje en covariate (uafhængig variabel) i stedet for 1 (intercept-only model).
Call: survfit(formula = Surv(dur, status) ~ udd, data = df)
udd=0
time n.risk n.event survival std.err lower 95% CI upper 95% CI
2 472 1 0.9979 0.00212 0.9937 1.000
3 469 2 0.9936 0.00367 0.9865 1.000
4 467 2 0.9894 0.00473 0.9801 0.999
5 465 2 0.9851 0.00558 0.9742 0.996
6 463 1 0.9830 0.00596 0.9714 0.995
8 462 2 0.9787 0.00665 0.9658 0.992
11 460 2 0.9745 0.00727 0.9603 0.989
12 458 7 0.9596 0.00908 0.9419 0.978
13 451 1 0.9575 0.00931 0.9394 0.976
14 450 3 0.9511 0.00995 0.9318 0.971
16 447 6 0.9383 0.01110 0.9168 0.960
17 441 1 0.9362 0.01127 0.9143 0.959
19 439 3 0.9298 0.01179 0.9070 0.953
20 436 6 0.9170 0.01273 0.8924 0.942
22 430 1 0.9149 0.01288 0.8900 0.940
23 429 3 0.9085 0.01331 0.8827 0.935
24 426 6 0.8957 0.01411 0.8684 0.924
25 420 2 0.8914 0.01436 0.8637 0.920
26 418 11 0.8679 0.01563 0.8378 0.899
27 407 2 0.8637 0.01584 0.8332 0.895
28 405 8 0.8466 0.01663 0.8146 0.880
29 397 2 0.8423 0.01682 0.8100 0.876
30 395 6 0.8296 0.01736 0.7962 0.864
31 389 3 0.8232 0.01761 0.7893 0.858
32 386 8 0.8061 0.01825 0.7711 0.843
33 378 1 0.8040 0.01833 0.7688 0.841
34 377 1 0.8018 0.01840 0.7666 0.839
35 376 2 0.7976 0.01855 0.7620 0.835
36 374 5 0.7869 0.01890 0.7507 0.825
37 369 2 0.7826 0.01904 0.7462 0.821
38 367 4 0.7741 0.01931 0.7372 0.813
39 363 1 0.7720 0.01937 0.7349 0.811
40 362 10 0.7506 0.01997 0.7125 0.791
41 352 3 0.7443 0.02014 0.7058 0.785
42 349 5 0.7336 0.02041 0.6947 0.775
43 344 5 0.7229 0.02066 0.6835 0.765
44 339 5 0.7123 0.02090 0.6725 0.754
45 334 3 0.7059 0.02104 0.6658 0.748
47 331 3 0.6995 0.02117 0.6592 0.742
48 328 13 0.6717 0.02168 0.6306 0.716
49 315 1 0.6696 0.02172 0.6284 0.714
50 314 6 0.6568 0.02192 0.6152 0.701
52 308 28 0.5971 0.02265 0.5543 0.643
53 280 2 0.5928 0.02269 0.5500 0.639
54 278 2 0.5886 0.02272 0.5457 0.635
56 276 3 0.5822 0.02277 0.5392 0.629
57 273 1 0.5800 0.02279 0.5371 0.626
58 272 3 0.5736 0.02284 0.5306 0.620
59 269 1 0.5715 0.02285 0.5284 0.618
60 268 7 0.5566 0.02294 0.5134 0.603
61 261 2 0.5523 0.02296 0.5091 0.599
62 259 3 0.5459 0.02299 0.5027 0.593
64 256 5 0.5353 0.02303 0.4920 0.582
65 251 1 0.5331 0.02304 0.4898 0.580
66 250 1 0.5310 0.02304 0.4877 0.578
68 249 3 0.5246 0.02306 0.4813 0.572
70 246 1 0.5225 0.02306 0.4792 0.570
72 245 6 0.5097 0.02308 0.4664 0.557
73 239 1 0.5075 0.02309 0.4643 0.555
74 238 1 0.5054 0.02309 0.4621 0.553
75 237 1 0.5033 0.02309 0.4600 0.551
77 236 1 0.5011 0.02309 0.4579 0.549
78 235 11 0.4777 0.02306 0.4346 0.525
79 224 2 0.4734 0.02306 0.4303 0.521
80 222 5 0.4628 0.02302 0.4198 0.510
81 217 1 0.4606 0.02302 0.4177 0.508
82 216 5 0.4500 0.02297 0.4071 0.497
83 211 1 0.4478 0.02296 0.4050 0.495
84 210 5 0.4372 0.02291 0.3945 0.484
85 205 1 0.4350 0.02289 0.3924 0.482
87 204 1 0.4329 0.02288 0.3903 0.480
88 203 3 0.4265 0.02284 0.3840 0.474
90 200 2 0.4222 0.02281 0.3798 0.469
91 198 1 0.4201 0.02279 0.3777 0.467
92 197 2 0.4158 0.02276 0.3735 0.463
94 195 1 0.4137 0.02274 0.3715 0.461
95 194 1 0.4116 0.02272 0.3694 0.459
96 193 4 0.4030 0.02265 0.3610 0.450
98 189 2 0.3988 0.02261 0.3568 0.446
99 187 1 0.3966 0.02259 0.3548 0.443
100 186 4 0.3881 0.02250 0.3464 0.435
102 182 3 0.3817 0.02243 0.3402 0.428
104 179 10 0.3604 0.02217 0.3195 0.407
105 169 1 0.3583 0.02214 0.3174 0.404
106 168 2 0.3540 0.02208 0.3133 0.400
108 165 2 0.3497 0.02202 0.3091 0.396
109 163 1 0.3476 0.02199 0.3070 0.393
110 162 4 0.3390 0.02186 0.2987 0.385
112 158 3 0.3325 0.02176 0.2925 0.378
114 155 1 0.3304 0.02173 0.2904 0.376
116 154 2 0.3261 0.02166 0.2863 0.371
118 152 1 0.3240 0.02162 0.2842 0.369
120 151 2 0.3197 0.02154 0.2801 0.365
122 149 1 0.3175 0.02151 0.2781 0.363
124 148 1 0.3154 0.02147 0.2760 0.360
129 147 1 0.3132 0.02143 0.2739 0.358
130 146 5 0.3025 0.02122 0.2636 0.347
131 140 1 0.3003 0.02118 0.2616 0.345
132 139 3 0.2939 0.02105 0.2554 0.338
133 136 1 0.2917 0.02101 0.2533 0.336
134 135 1 0.2895 0.02096 0.2512 0.334
135 134 2 0.2852 0.02087 0.2471 0.329
142 131 4 0.2765 0.02069 0.2388 0.320
144 127 3 0.2700 0.02054 0.2326 0.313
145 123 1 0.2678 0.02049 0.2305 0.311
146 122 2 0.2634 0.02039 0.2263 0.307
148 120 4 0.2546 0.02017 0.2180 0.297
150 116 2 0.2502 0.02006 0.2138 0.293
152 113 4 0.2414 0.01984 0.2055 0.284
154 108 1 0.2391 0.01978 0.2034 0.281
156 105 10 0.2164 0.01916 0.1819 0.257
158 83 2 0.2111 0.01905 0.1769 0.252
160 77 3 0.2029 0.01889 0.1691 0.244
161 67 1 0.1999 0.01885 0.1662 0.240
162 66 2 0.1938 0.01876 0.1603 0.234
163 63 1 0.1908 0.01871 0.1574 0.231
165 57 1 0.1874 0.01868 0.1542 0.228
166 56 1 0.1841 0.01864 0.1509 0.224
168 52 1 0.1805 0.01862 0.1475 0.221
182 46 3 0.1688 0.01860 0.1360 0.209
183 37 1 0.1642 0.01865 0.1314 0.205
194 24 1 0.1573 0.01909 0.1241 0.200
208 22 2 0.1430 0.01985 0.1090 0.188
210 16 1 0.1341 0.02053 0.0993 0.181
220 9 1 0.1192 0.02303 0.0816 0.174
234 5 1 0.0954 0.02818 0.0534 0.170
352 2 1 0.0477 0.03654 0.0106 0.214
364 1 1 0.0000 NaN NA NA
udd=1
time n.risk n.event survival std.err lower 95% CI upper 95% CI
7 77 1 0.987 0.0129 0.962 1.000
28 76 1 0.974 0.0181 0.939 1.000
30 75 1 0.961 0.0221 0.919 1.000
31 74 1 0.948 0.0253 0.900 0.999
34 73 2 0.922 0.0305 0.864 0.984
35 71 1 0.909 0.0328 0.847 0.976
36 70 1 0.896 0.0348 0.830 0.967
48 69 1 0.883 0.0366 0.814 0.958
51 68 1 0.870 0.0383 0.798 0.949
52 67 1 0.857 0.0399 0.782 0.939
54 66 1 0.844 0.0413 0.767 0.929
56 65 1 0.831 0.0427 0.752 0.919
57 64 1 0.818 0.0440 0.736 0.909
68 63 2 0.792 0.0462 0.707 0.888
69 61 1 0.779 0.0473 0.692 0.878
76 60 1 0.766 0.0482 0.677 0.867
78 59 1 0.753 0.0491 0.663 0.856
86 57 1 0.740 0.0500 0.648 0.845
96 56 2 0.714 0.0516 0.619 0.822
100 54 1 0.700 0.0523 0.605 0.811
103 53 1 0.687 0.0530 0.591 0.799
104 52 4 0.634 0.0551 0.535 0.752
108 47 1 0.621 0.0555 0.521 0.740
109 46 1 0.607 0.0560 0.507 0.728
122 45 1 0.594 0.0563 0.493 0.715
128 43 1 0.580 0.0567 0.479 0.702
142 40 1 0.566 0.0571 0.464 0.689
143 39 1 0.551 0.0574 0.449 0.676
150 37 1 0.536 0.0578 0.434 0.662
151 36 1 0.521 0.0581 0.419 0.648
156 34 1 0.506 0.0583 0.404 0.634
182 25 1 0.486 0.0594 0.382 0.617
184 21 1 0.463 0.0609 0.357 0.599
187 13 1 0.427 0.0658 0.316 0.578
209 8 1 0.374 0.0762 0.250 0.557
212 7 1 0.320 0.0819 0.194 0.529
220 6 1 0.267 0.0839 0.144 0.494
332 1 1 0.000 NaN NA NA
Har forskellige grupper forskellig overlevelse?
Første skridt i undersøgelsen af hvorfor (i et kausalt henseende).
Vi tester:
\[H_{0}: S_{1}(t)=S_{2}(t)=\dots=S_{n}(t)\]
Med andre ord:
Nul-hypotese: de to overlevelsesfunktioner er ens. \(H_{0}: S_{1}=S_{2}\)
Alternativ hypotese: de to hypoteser er forskellige. \(H_{1}: S_{1} \neq S_{2}\)
Ved flere grupper er logikken den samme. Alle overlevelsesfunktioner er ens eller mindst to overlevelsesfunktioner er forskellige
Der er forskellige signifikanstests, men log-rank (Mantel-Haenszel) test er den mest anvendte.
Teknisk er en log-rank test en “large scale chi-squared test”, der estimerer en overordnet—global—sammenligning af to eller flere KM kurver.
Det gør vi ved at udvide de forrige tabeller—dem med tid, hændelse, risksæt, censureringer, osv.—til også at inkludere forventede antal hændelser. Forventede hændelser, \(e\), for hver gruppe, er givet ved:
\[ e_{1} = \left ( \frac{n_{1}}{n_{1}+n_{2}} \right ) \times \left ( d_{1}+d_{2} \right ) \] og
\[ e_{2} = \left ( \frac{n_{1}}{n_{1}+n_{2}} \right ) \times \left ( d_{1}+d_{2} \right ) \] Herfra kan vi udregne forskellen i observeret \(O\) og forventet \(E\), som:
\[ O_{i}-E_{i}=\sum \left ( d_{i} - e_{i} \right ) \]
hvor \(i\in1,2\).
Summen for de to grupper vil være det samme, bortset fra den ene har minus foran summen.
Log-rank statistikken for to grupper (1, 2) er slutteligt givet ved:
\[ \frac{\left ( O_{1} - E_{1}\right )^{2} }{\text{Var} \left ( O_{1} - E_{1} \right ) } \] eller
\[ \frac{\left ( O_{2} - E_{2}\right )^{2} }{\text{Var} \left ( O_{2} - E_{2} \right ) } \]
P-værdien udleder vi fra \(\chi^{2}\) fordelingen på klassisk vis. Dette lader vi selvfølgelig computeren gøre, og printes ofte i bunden eller toppen af outputtet, afhængigt af kode-sprog.
Alternative til log-rank er Wilcoxon, Tarone-Ware, Peto, og Flemington-Harrington tests.
Disse er alle variationer af log-rank testen med har forskellige vægte, der hhv. vægter tidlige og sene hændelser i forløbet. God praksis er at vælge en test a priori baseret på forventninger om dataen og fænoment. Disse forventninger kan trykprøves ved også at lave de alternative tests. Ofte vil de andre test have samme signifikans-niveau.
Call:
survdiff(formula = Surv(dur, status) ~ udd, data = df)
N Observed Expected (O-E)^2/E (O-E)^2/V
udd=0 474 390 348.2 5.01 26.3
udd=1 79 44 85.8 20.34 26.3
Chisq= 26.3 on 1 degrees of freedom, p= 0.0000003
single-episode file med fixed/time-invariant variable:single-episode file med fixed/time-invariant variable:PNR (individuel ID)aar (året for hændelsen)ALDER (alder på hændelsestidspunktet)KOEN (kønsvariabel – baseret på cpr)f_udd (højeste fuldførste uddannelse)region18 (bopælsregion da personen er 18 år)event (1=hændelse, 0=uændret tilstand)I skal undersøge tiden indtil første fødte barn. Dette eksempel fortsætter gennem forelæsningerne, således i kan bygge ovenpå hver øvelse. MEN, har i lyst til at arbejde med en anden problemstilling i finder interressant, eller evt. matcher med jeres aktuelle projektarbejde, er i også velkomne til det.
Diskuter hvorfor variablene er fixed/“tids-invariante”? Hvorfor ikke?
Konstruer en meningsfuld tids-variabel ud fra ALDER, således alle starter i populationen under risiko (tid = 1), når de fyldte 18 år.
3.1 Lav en KM estimation: beskriv og fortolk resultaterne.
3.2 Lav en meningsfuld og fyldestgørende visualisering.
region18) i tiden indtil første fødte barn.Call: survfit(formula = Surv(tid, event) ~ 1, data = First_child)
n events median 0.95LCL 0.95UCL
[1,] 3464 2322 12 12 12
Jeppe F. Larsen | 28. februar 2024